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Abstract: In this paper we study the AdS/CFT correspondence for = 4 SYM 
with gauge group U{N), compactified on in four dimensions using Monte-Carlo 
techniques. The simulation is based on a particular reduction of degrees of freedom 
to commuting matrices of constant fields, and in particular, we can write the wave 
functions of these degrees of freedom exactly. The square of the wave function is 
equivalent to a probability density for a Boltzman gas of interacting particles in six 
dimensions. From the simulation we can extract the density particle distribution for 
each wave function, and this distribution can be interpreted as a special geometric 
locus in the gravitational dual. Studying the wave functions associated to half-BPS 
giant gravitons, we are able to show that the matrix model can measure the Planck 
scale directly. We also show that the output of our simulation seems to match 
various theoretical expectations in the large limit and that it captures effects 
as statistical fluctuations of the Boltzman gas with the expected scaling. Our results 
suggest that this is a very promising approach to explore quantum corrections and 
effects in gravitational physics on AdS spaces. 
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1. Introduction 

Many approaches to quantum gravity propose to replace the continuous geometry of 
spacetime by a discrete approximation whose underlying discreteness is concentrated 
at the Planck scale, the distance scale at which quantum effects are important. It is 
assumed that this discrete theory is the correct description of gravitational physics 
at the Planck scale (recent examples of such ideas in concrete settings amenable 
to computer simulations can be found in f^). The main technical problems with 
these approaches is to show that in the continuum limit (if it exists) one reproduces 
Einstein's theory of general relativity. This continuum limit would be related to 
quantization of gravity on large manifolds (as measured in Planck units), and classical 
geometry would be emergent in the sense that the continuum geometry is not part 
of the discrete problem, but that it only appears in a suitable limit. This will only 
be possible if one has a consistent quantum mechanical evolution (thereby solving 
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the problem of time), as well as a semiclassical state of the theory that permits 
making the comparison to a semiclassical or classical approach based on perturbative 
techniques ^. 

Once one would obtain such a theory, it might be very hard to solve it in general 
and it is desirable that the discrete theory could be simulated on a computer. The 
main objective would be to try to analyze typical configurations and to ask what 
new phenomena can be understood in this way. Clearly, if one wants to simulate an 
infinite space like Minkowski space or AdS geometries, one would need an infinitude 
of points and then such type of computational approach would need to be regulated 
somehow. 

There are other approaches to a definition of quantum gravity that do not be- 
gin by assuming that spacetime is there as a fundamental object either. Within 
string theory, a very concrete proposal for a formulation of quantum gravity on AdS 
geometries that implicitly solves the problem of time is given by the AdS/CFT cor- 
respondence [0. The simplest example of the AdS/CFT correspondence states that 
type IIB string theory on an asymptotically AdS^ x spacetime with units of 
flux is exactly equivalent as a quantum mechanical system to the maximally super- 
symmetric SYM theory with gauge group SU{N) in four dimensions, compactified 
on the conformal boundary of (global) AdS^ x . This conformal boundary is 
X M. Thus, the quantum field theory is compactified on a round three sphere, and 
one has time translation invariance associated to the M direction. Since the SYM 
theory has a consistent quantum mechanical evolution and Hamiltonian, the problem 
of gravitational physics is to understand how the geometry of the ten dimensional 
spacetime is encoded in the quantum system associated to the field theory. This can 
be rephrased as the following question: how does metric geometry emerge? 

In particular, the SYM theory should in principle contain all of the (consistent) 
solutions of Einstein's equations with the required boundary conditions, so long 
as their radius of curvature is large in Planck units everywhere. These classical 
gravity solutions should be associated to some semi-classical state in the gravitational 
formulation, and for each of these one would have a dictionary to a particular state 
of the SYM theory (this is, a non-trivial wave function of the SYM theory). The 
dynamical evolution of these CFT states should reproduce in some way the dynamical 
evolution of the gravitational solution. 

In this paper we want to show that many^ aspects of this question can be analyzed 
numerically in the CFT, whereby one would recover a lot of information of various ten 
dimensional metrics from doing a computer simulation with finitely many degrees of 
freedom. The simulation we propose and execute in this particular paper is a Monte- 

^To show that one can excite gravitational waves, or other quanta of the theory, one would need 
a family of such solutions, which are to be considered as small perturbations of the original solution. 

^The word many in this phrase can be replaced by some or few, according to the detailed 
questions one wants to address. 
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Carlo algorithm that explores wave functions of an associated problem of N identical 
particles in six dimensions. The relations of this particular system to the strong 
coupling limit of the CFT were explained and justified in and we will revisit them 
later on in the introduction. Given the wave functions, the distribution of particles 
in six dimensions in the thermodynamic limit (moderately large N) is concentrated 
on a five dimensional submanifold of M^, and this submanifold is associated to a 
five-dimensional slice of the geometry of the ten dimensional spacetime. 

In this paper we will be testing aspects of this idea with the simplest possible 
configurations, and we will show that there is an opportunity to study aspects of 
quantum gravity that are not accessible otherwise with these numerical methods. 
The main goal of this paper is to explain how such an identification is possible and 
to give numerical evidence supporting such an identification. 

From a more general perspective, one of the biggest problems we have to face 
when trying to match gravity solutions and SYM states is that the SYM states have 
to be evaluated at strong coupling, and very little is known about the behavior of 
SYM at strong coupling from first principles. In a similar vein, very little is known 
about the "most general solution of Einstein's equations" with the requisite boundary 
conditions. We need a simpler problem, just because our calculational abilities are 
limited. 

Thus, it seems reasonable to restrict the problem further and exploit the super- 
symmetry of the supergravity and CFT problems to get to a more manageable set 
of states to compare. In this sense, we can choose to find states that preserve many 
supersymmetries and symmetries. With these restrictions it is possible to make some 
headway into the problem of matching both sides of the duality. This is the main 
simplification of the problem that lets us reduce the degrees of freedom sufficiently 
to make a detailed comparison possible. In principle one could expect that other 
simplifications that use symmetries can result in a similar reductions of the degrees 
of freedom. 

It seems reasonable to choose states that preserve 1/8 of the supersymmetries. 
The reason to choose such states, is that such a class of states can be associated to 
the chiral ring of the SYM via the operator state correspondence, and the chiral ring 
has many non-renormalization theorems that permit a sensible comparison. This is 
also done to avoid the potential problem of having to evaluate a very complicated set 
of quantum corrections. The other advantage of choosing such states, is that they 
are guaranteed to be associated to supergravity solutions of string theory, because 
the only one-particle states saturating the BPS bound in the gravitational theory 
are those associated to massless particles [^, and the only massless particles lie 
in the supergravity multiplet. Thus, in this case we do not have to test if we have a 
state made of heavy strings before comparing to geometry: one is guaranteed to be 
testing solutions of the type II supergravity theory in ten dimensions. 

Our analysis of BPS states begins with the field theory (CFT) setup. It is 
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convenient to work in the canonical formalism of the field theory dynamics, for the 
field theory compactified on x M. The SYM theory multiplet contains as part of 
its spectrum 6 scalar fields (f)°'{x,t) in the adjoint of U{N), where x is the position 
of the field in (the dynamics of the U{1) degrees of freedom are decoupled and 
can be ignored). To solve the theory in the free field limit, one would decompose the 
fields into spherical harmonics of the sphere S"^, in a normal mode decomposition 

rix,t) = J2^lmmirnix) (1.1) 

Here we are abusing notation as if the 5*^ where a two-sphere. ^ The important point 
is that we have now a sum over a discrete set of oscillators, and it is straightforward 
to compute their energies, wi^rn ~ / + 1. 

A straightforward analysis of the chiral ring states shows that the only oscillators 
that are turned on near the free field limit are those with / = [^], what one 
would call the s-wave of the normal mode expansion. These are just constant field 
configurations. 

Because of these facts, and because the states we want to explore are super- 
symmetric, one can assume that even at strong coupling only such constant field 
configurations are relevant, and that supersymmetry will guarantee that a semiclas- 
sical approach is essentially exact [H. 

The advantage of this formulation for a partial solution via a BPS problem, is 
that it reduces the problem of infinitely many degrees of freedom of the field theory 
(all the spherical harmonics of the sphere) to finitely many degrees of freedom. 

For these degrees of freedom, we can write an effective Hamiltonian, which is the 
dimensional reduction of the field theory on a sphere to the constant modes. This 
hamiltonian is given by 

^e// = tr (X^^(AXO^ + ^(XT+ E ^[X^X^][X^X^]] (1.2) 

\i=l a,b=l J 

where the X* are six N x N hermitian matrices, and we have gauge transformations 
X' — s> U^^X^U where the U are unitary matrices. The factors of Svr^ in the inter- 
actions arise from choosing canonical normalization of the matrix fields, once the 
volume of the three sphere is taken into account. A further reduction in degrees of 
freedom is obtained from going towards strong 't Hooft coupling, i.e. Qym^ large. 
This is because we usually expect that for random matrices, the eigenvalues are of 
order \fN . Thus the term that is multiplied by Qym g^ts an extra factor of X^, 
which tells us that at large X, the potential term dominates the dynamics. 

^In practice, the has an SU{2) x SU{2) isometry group, and the normal modes of a scalar 
appear in the (n/2,n/2) representations of this group. Thus, the principal quantum number n 
classifies the states with the same energy. In our abused notation, I ~ n and m includes the spin z 
quantum numbers of both SU (2) groups. 
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The idea is that to determine the ground state, as well as the low lying excitations 
of the system we should look at the matrix configurations that minimize the quartic 
term, and localize on those configurations. These configurations require the matrices 
X* to commute. There is still a family of matrix configurations that have these 
properties. 

This is a very non-trivial statement, but there is strong evidence that this 
statement is correct, by comparing the energies of semi-classical string states [[7| 
in AdS^ X with the results of such an approach (see also p), where the string 
states are constructed by turning on the off-diagonal modes of the X fields in a 
perturbative expansion. These also confirm other perturbative resumations of pla- 
nar diagrams |jll| and input from integrability structures by doing a calculation 
directly at strong coupling. 

Being more careful about supersymmetry, one can notice that these configura- 
tions of commuting matrices solve the F and D-term constraints of the field theory 
on fiat space and are associated to the moduli space of vacua of the Af = 4 SYM. 
Having commuting matrices is required in order to preserve supersymmetry, but it is 
not the only condition that one needs. One also needs a holomorphic wave function 
[0]. Via a gauge transformation, we can choose all the X* to be diagonal simultane- 
ously for these configurations, and then the dynamics of supersymmetric states in the 
chiral ring is a particular quantization of the moduli space of vacua of A/" = 4 SYM. 
This moduh space is described by N (identical) particles on C'^ (or M^). Because 
we also have a kinetic term and a quadratic term, we are going to obtain particular 
normalizable wave functions on this moduli space, and there will also be a preferred 
origin where X = 0. The coordinates of these particles are the different (correlated) 
eigenvalues of the X* (this is in the same spirit as Matrix theory [|l2l). Thus we can 
associate a six vector position per eigenvalue, as follows 



X 
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(1.3) 



The quantum mechanics obtained in this way is almost free, as we have removed 
the explicit interactions between the matrices by requiring them to commute. If we 
diagonalize the matrices, we get an extra measure term /i^ from going to diagonal 
variables This has been generalized to some orbifolds [|T3[, and suggests that 
one can generalize most of our techniques to the orbifold setup without too much 
trouble. 

The problem of the study of the full CFT is reduced in this way to a study of a 
dynamical problem of N particles in six dimensions with a given Hamiltonian. All we 
have to do now is solve the dynamics and find a large class of wave functions ip that 
solve the (time independent) Schrodinger equation. Given any such wave function 
ip, we then take a modified wave function ip = fiip. The reason why it is convenient 
to absorb fi in the definition of ip, is that fi^ will be the measure that is required for 
orthogonality of different eigenvalues of the Hamiltonian H. 
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From here is a probability density in the usual sense, and it describes a gas 
of correlated particles on C^. It is exactly the distributions of these particles in 
six dimensions that are interesting for our purposes. The reason for this, is that one 
is supposed to recover gravity in the classical limit as we send N ^ oo, and one 
expects the gravitational description to be some form of collective coordinates or 
thermodynamic description of the field theory problem. The simplest such collective 
coordinate is the density distribution of particles in this six dimensional space of a 
"single particle" location. 

It turns out that one can find a lot of approximate solutions of the Schrodinger 
equation for the Hamiltonian obtained by this procedure, and that these are all 
normalizable wave functions. The type of wave function that one obtains is similar to 
a thermal gas of particles with a confining potential and with logarithmic repulsions 
in six dimensions, which are induced by taking the logarithm of //^. This type 
of problem is especially suited to Monte-Carlo simulations. This is very easy to 
implement on a computer. In this paper we will explain the algorithm we use to 
explore these configurations, as well as some of the improvements required to take 
large in a systematic way. We will also check our numerical results with the 
known theoretical results for this kind of thermodynamic gas in the large limit 
(unfortunately very little is known analytically). 

The connection of this (commuting) matrix quantum mechanics model to gravity 
requires us to use the ideas of the AdS/CFT correspondence in detail in order to 
extract gravitational information from these configurations. We will explain our 
proposal for what kind of gravitational information we can expect to extract from 
these distributions of particles. By comparing with other known behavior of exact 
BPS solutions in gravity [l^, we will be able to see that different topologies of 
spacetime are available to us. 

Our paper is organized as follows. In section ^ we explore the Hamiltonians 
that we are required to analyze, as well as the wave functions that we will study. 
We will also describe and analyze some properties of the these wave functions in 
the thermodynamic limit (large N limit). Some more results that are relevant for 
this section but are not required for the main argument are in appendix ^ This 
section can be skipped on a first reading. In section ^ we describe the chiral ring 
and it's relation to gravitational solutions. We discuss in detail the half-BPS states 
and we review briefly the work of Lin, Lunin and Maldacena \TM. We use analogy to 
establish a relation between the geometry of the particle distributions obtained from 
the wave function calculations and the geometry of the dual gravitational solution. In 
section ^ we describe in detail our Monte Carlo algorithm as well as how we propose 
to handle the large A^ limit. Our numerical results are compared to analytical large 
A^ theoretical results in section |^. We show that the numerical calculation agrees 
well with the theoretical expectation for the matrix models. We also study some 
configurations that are dual to LLM metrics, and that are associated to different 
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topologies of spacetime. We see a clear topology change. Moreover, we are able to 
show that the matrix model we are studying can be used to infer the size of the 
geometry in Planck units. This is, we are seeing evidence for effects that measure 
the Planck length in the simulation. On the other hand, some of the features of the 
simulation do not seem to agree with the LLM intuition. We end the paper with our 
conclusions and the prospects for future simulations along the lines we have outlined. 



2. Relevant exact results for quantum matrix models of com- 
muting matrices 

As described previously, we have reduced the problem of the full CFT dynamics of 
Af = 4 SYM to a simple matrix quantum mechanics, and we further reduced that 
problem to a gauged quantum mechanical model of commuting matrices. 

Because we can use gauge transformations to diagonalize matrices, the eigenval- 
ues of the matrices are the only relevant degrees of freedom. However, permutations 
of the eigenvalues can be obtained by choosing to diagonalize the matrices via a 
different unitary operator. As such, the permutations of the eigenvalues are a gauge 
symmetry, which can be implemented by requiring that the wave function is invariant 
under such permutations. This gives us a problem of bosons in a non-trivial setup. 

Our purpose in this section is to explore the following Hamiltonian for N identical 
(bosonic) particles in d spatial dimensions 

^ = E-^^^-'"'^^ + ^l^^l' (2.1) 



where 

,2 



/i 

i<j 



Y[\xi-x,\^. (2.2) 



in the case where the dimension d is even (and eventually it will be strictly greater 
than two). This Hamiltonian was derived in |^. We refer the reader to that paper 
for details. 



We should notice that if fx were replaced by one in the Hamiltonian above 
would look just like that of N particles in a d-dimensional harmonic oscillator, this 
is 



It is also clear that /i can only be important where | small. It seems 

reasonable to try to solve the Hamiltonian with the Gaussian wave function that 
would solve the simple harmonic oscillator, plus some corrections where Xj — Xj is 
small. 



-7- 



If we try 

V'o~exp(-^f2/2) (2.4) 

i 

we can ask if ipQ is actually an eigenstate of the above Hamiltonian. 
Direct manipulation shows that 

ViV'o = Xiipo (2.5) 

Similarly we can calculate 

V, ■ (/i'xO^o = (x^) W^o + (Vi ■ Xi)/i2V;o + iJoXi ■ Vifi^ (2.6) 
from which we can calculate 

i 

easily 

We notice the following facts. The first term behaves like xfipo and cancels the 
quadratic potential piece in the Hamiltonian. The second term is a constant times 
ipQ. (This constant is Nd). The third term is 

^oJ^^i-V./i^ (2.7) 

i 

We recognize that Xj ■ Vj is an Euler vector for dilatations where Xj Xxi 
simultaneoulsy. Clearly /x^ is an eigenvalue of this operator with eigenvalue A^(A^— 1), 
because /i^ is a homogeneous function of the Xi of this degree. 
As such we obtain that 

Htfjo = ^0^0 (2.8) 

where Eq = N{N — l)/2 + Nd/2. It should be pointed out that in the case d = 1, 
this is exactly the energy of A^^ harmonic oscillators, just as expected, because this 
is a single matrix model quantum mechanics. This was solved in (the review 



article |T^ is also very useful to get some intuition of these setups). 

It is convenient at this point to define a new wave function by a similarity 
transformation -ijj = flip. The measure of ip that makes the Hamiltonian self-adjoint 
is given exactly by /i^. This is, the probability for the particles to be located at 
various positions (within a fixed region R about a particular set of positions) is 



fi^ilj*iP ~ / iP*iIj (2.9) 

R Jr 

Let us change the d real coordinates of Xj, x^ . . . by the following complex 
combinations 

z'. = xf-^ + ixf (2.10) 
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by choosing a complex structure for the d dimensional space. 

Now, let us consider a homogeneous polynomial of degree m in the Zj coordinates, 
and let us call this polynomial P. Furthermore, we need P to be symmetric when 
we exchange the particles in C^. With this additional information, we construct the 
wave function 

i;p = ^oP{z]) (2.11) 

The conjecture stated in |^ is that all of these wave functions are the exact 
wave functions of the dual states to the chiral ring, and that the energy of the 
configuration is given exactly by the degree of the polynomial. In the appendix ^ 
we study these wave functions with the naive Hamiltonian written above where we 
have broken manifest supersymmetry by our choices of which degrees of freedom we 
are keeping. One can state that these wave functions are very good approximations 
to exact solutions of the non-supersymmetric Hamiltonian above. So if one were to 
do the same analysis being careful about supersymmetry, it is likely that one might 
be able to prove exactness of the wave functions. It would be desirable to have such 
a formulation. For our purposes, all we need is the knowledge of the wave functions 
themselves. 

Their R charge (quantum numbers under rotations of the X fields) of the wave 
function is also the degree. 

Notice also that when we calculate an overlap of two holomorphic wave functions 
with different degrees in the form 

j \tPo\'Pm{znPm'{z') (2.12) 

we necessarily get zero. This is because the measure and is invariant under phase 
rotations — > exp{i6)z'' if we do it on all the z simultaneously, while Pm{z) — > 
eyi\i{im6)Pm{z). This is related to a particular diagonal 50(2) charge inside the 
5*0(6) R-charge of A/" = 4 SYM theory, which defines an A/" = 2 polarization. 

Compared to the vacuum, these states have energy m, provided we redefine the 
zero of energy so that i?o = 0. These states also have i?-charge m. This is, the 
wave functions we have constructed have essentially the same energy and i?-charge 
as expected for a BPS multiplet. Because A/" = 4 SYM is conformally invariant, one 
can also use the operator-state correspondence. This correspondence states that for 
every local gauge invariant operator of the theory in (Euclidean) flat space one can 
find a quantum state \0) of the field theory compactified on . The correspondence 
further states that the conformal dimension of O is the same thing as the energy 
of \0) with respect to the ground state of the theory on 5^. The set of states 
we have obtained would have the same i?-charge and conformal dimension (energy) 
as elements of the chiral ring. Structurally, the chiral ring is also built by gauge 
invariant homolorphic operators, where the F term constraints are imposed. These 
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are holomorphic functions on the moduh space of vacua of the theory, and they would 
correspond exactly to symmetric holomorphic polynomials like P. 

For the rest of the paper we will concentrate on studying the holomorphic wave 
functions described above and we will treat them as exact wave functions. It would 
be interesting if this could be proved exactly by using supersymmetry arguments. 

2.1 Thermodynamic behavior of the A^-particIe wave functions 

We now have a list of wave functions to analyze. They are all built by multiplying 
the ground state wave function %pQ by a symmetric polynomial of the variables 2;*. 
We now want to find out what type of geometry these wave functions are associated 
to. 

To begin, we want to study the ground state itself. We find that the square of 
the wave function, which has a probabilistic interpretation, takes the following form 



Notice the similarity between this probability function and a Boltzman distri- 



^j^i — X]j<j^og(l^i ~ ^iP; where the Xj are the positions of the particles. This is, 
we notice that we have a gas of particles confined by a harmonic oscillator well, and 
that also have repulsive logarithmic interactions. If ci = 2 this is a Coulomb gas of 
particles (a two dimensional plasma) in a potential well. For higher dimensions this 
is a different problem. 

If we are interested in a thermodynamic limit where N 00 (meaning N is 
taken to be very large), then we can hope that the gas will settle to a preferred ther- 
modynamic configuration that will maximize the probability distribution, and that 
is well described by a density of particles in d dimensions with some characteristic 
thermal fluctuations. 

We would approach this thermodynamic limit replacing sums by integrals, and 
introducing a density of particles p{x). In this way the energy is given by 




(2.13) 



i<j 




where /3 



1 and H 




(2.14) 



subject to 





and to p{x) > 0. As is typical in thermal problems, we would first find the saddle 
point of H that minimizes the energy. This was done in detail in The main 
observation is that we need to solve the following integral equation 




(2.16) 
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in the region where p has smooth support, and where C is a lagrange multipher 
enforcing the constraint. If d is even, and greater than two, then 



iylf'\og{\x-y\)^5\x-y) 



(2.17) 



so that under the assumption of smooth support, taking derivatives with respect to x 
and integrating over y commute. If one uses this assumption, one applies {{V xYY^"^ 
on both sides of the equation. After this procedure one obtains the following equation 



This contradicts that J p = N. So one must conclude that the particle distribution 
p has singular support in the thermodynamic saddle point limit. For the case above, 
it was found that due to spherical symmetry, one expects that the distribution is 
spherically symmetric and singular p{x) ~ S{\x\ — r). The saddle of this ansatz 
occurs exactly when 



swer. As we will see in section ^, our simulation proves that this is correct in the 
thermodynamic limit. 

For the case of A/" = 4 SYM, we need the special case d = 6. In this case there 
are three (matrices of) complex variables that we need. It is customary to call them 
X, Y, Z. We have already used X" as real variables. The notation we will follow is 
that X without an index represents X^ + iX"^. There are other simple symmetric 
polynomials P{X, Y, Z) one might consider other than one. 

For example, take a single trace polynomial of Z, 



the wave function ijjoPn is an allowed wave function. So is ipi = ?/'o(tr(Z"'))^ = 
ipoPn and ^2 = ipo^^iZ'^"') = "^o-Pan- These two states are distinct functions of the 
eigenvalues of Z with the same energy and R charge, so they represent different states. 
What can we say about these two wave-functions? Here we need some more intuition. 
Since we are interested in comparing these wave functions with a gravitational dual 
configuration of AdS, we will get our intuition of what these objects represent from 
the expected dictionary of the AdS/CFT setup. 

In the AdS/CFT dictionary established by Witten 0, each trace counts as a 
single graviton, so one expects that at large N, for fixed n, the two states ipi and 
ijj2 (properly normalized) are approximately orthogonal. This is, {ipi\'!^2) ~ O^N"^). 
The number of traces counts approximately the number of gravitons, and one can 
check that the traces do have an approximate oscillator algebra for a single matrix 







(2.20) 
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model, where ajjct) ~ Pn{Z)iJja. Here we have to assume this property, as we don't 
know how to calculate the norm of the corresponding states analytically. 

Single graviton states on a given geometry do not correspond to classical states of 
geometry. To get geometrical states, one would naturally expect that these are given 
by some type of coherent state. This is, our first guess for interesting geometric wave 
functions is to take an expression of the form exp(^^i„a''^)|0) with finitely many tn 
different from zero. 

However, we see that we have a problem with a naive extrapolation of Witten's 
result for this type of wave function. The reason is that 

eMj2^nPniZ))iPo (2.21) 

n 

is not a normalizable state unless tn = for all n > 3, because the trace dominates 
over ipo for very large values of the eigenvalues of Z, z. Thus we can not do quantum 
mechanics with such a state. What we need to fix this is good behavior at infinity, 
so that the wave function is integrable. 

Let us define f{x) — Y^^^tnX"'. With this convention, Ylin^nPn{Z) — tr/(Z). 
To cure the bad behavior at infinity we can require that / have better behavior at 
infinity than any term in the series expansion of /. However, / is a complex analytic 
function, so to have such a property, / has to be multivalued. 

The simplest behavior is for f{x) to behave logarithmically at infinity. This is 
just thinking of a particular polynomial behavior of '?/'/'?/'o, and we know that such 
wave functions are normalizable. The function / will then have a branch cut in the 
complex plane. However can still be single valued, as it depends on exp(i/). Let 
us use instead 

f^m\og{g{x)) (2.22) 

where g is an arbitrary polynomial of x, and m is a parameter that tunes the strength 
of the perturbation. For single valued wave functions, it should be an integer. For x 
small, if 5^(0) ^ 0, we can expand \og{g{x)) in Taylor series, and we can approximate 
any polynomial / to arbitrary accuracy. Such a wave function would be given by 

Viodet(^(Z))- (2.23) 

We expect that these wave functions are the ones that are important for geome- 
try. Let us now consider a case where g{Qi) ^ 0, and where the first zero of g happens 
for X very large (much larger than \/iV, the typical radius of the sphere) . 

One would find then that the numerical value of / is small for x ~ V^, and that 
/ is holomorphic in this region. As such, we can think of / as a small perturbation of 
the confining potential, and that the saddle point of p will react by a small change. 

We can again go to the thermodynamic limit, and we find that we now need to 
satisfy the following integral equation 

x^ + C- f{z) - m = J d'^ypiy) log(|f - (2.24) 
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where z = + ix^. Thus V^(/(z) + f{z)) = 0, because / is holomorphic. Again, 
assuming p is smooth leads to a contradiction. This leads us to p being given by a 
singular distribution as well. The hard problem is to determine the distribution p 
given /. For situations where / is considered to be small, we expect that there is no 
change of topology, but that the geometric profile of the sphere gets deformed. 

A more general analysis shows that 

^odet{g{X,Y,Z)r (2.25) 

will behave the same way for g a polynomial in three complex variables, where the 
closest zero of g to the origin is at distances much larger than -\/iV- The factor of m is 
an integer controlling the strength of the perturbation of the ground state geometry. 

One can also consider situations where the zeroes of g are more generic and let 
us say that these zeroes intersect the ^. One finds that the particles are repelled 
from the zeroes of g, because the wave functions vanish there. Thus the support of p 
does not intersect the zeroes of g, and the same arguments apply: / is holomorphic 
(and multivalued) in the region of interest, but / + / is a real single valued function 
whose Laplacian vanishes. We expect in general that the description of all of these 
states corresponds to five- dimensional submaniolfds of of different topologies. 

It is interesting to ask the following qualitative questions: what are the allowed 
topologies of these submanifolds (let us say for fixed degree of g)7 Do these mani- 
folds have boundaries, or lower dimensional pieces? Are there universal features of 
topology transitions as we vary g with fixed degree? 

More quantitatively, we would like to know how the geometry of the submanifold 
correlates with the exact details of the wavefunction. Particularly, how the density of 
particles in the submanifold is related to other geometric properties of the embedding 
and if there is interesting scaling near topology transitions. One would also like to 
understand how finite effects blur the topology transitions. 

Our purpose in this paper is to show that these quantitative and qualitative 
questions about the particle distributions can be addressed numerically by simulating 
the wave functions using a Monte Carlo algorithm. 

Once these questions are understood (numerically) in the matrix model, one 
would like to find the gravitational dual description of these questions. Our proposal 
is that the geometric features discussed above correlate directly with the gravitational 
dual description and provide answers to quantum gravity questions for which there 
is no analytic understanding of finite A^ effects. We will explain this correlation in 
the following section. 

^Determinant operators are also related to giant gravitons in the AdS/CFT setup Thus 
one can identify the integer m above with the number of Z3-branes in this case 
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3. AdS duals to chiral primaries 



Chiral primary operators by definition are BPS objects. Their main property is that 
their conformal dimension A is equal to their R-charge, J. These operators are re- 
lated to BPS states for the field theory compactified on S^. The superconformal 
algebra guarantees that J and H commute. One can show that the twisted Hamil- 
tonian H = H — J is positive semi-definite in the free field limit and that only BPS 
states are annihilated by H. The classical solutions of if = are given exactly by 
constant configurations of the fields on and they are identical to the moduli space 
of vacua of the AT = 4 SYM theory [§. 

On the gravitational side, one thinks of J and A as particular isometries of 
AdS^ X [0] . For a point particle moving in this space A is the energy and J is the 
angular momentum of the particle motion. Setting A = J tells us that the energy is 
equal to the momentum. This is, the particle is massless in ten dimensions. For type 
IIB string theory, all massless particles belong to the supergravity multiplet under 
supersymmetry transformations. Thus all BPS configurations with this amount of 
supersymmetry and with small quantum numbers should correspond to a gravita- 
tional excitation of the AdS geometry. This is, one should be able to understand all 
these supersymmetric configurations in terms of supergravity. Because the different 
quanta respect the same supersymmetry, in principle it is possible that the gravi- 
tational attraction between two quanta is compensated exactly by exchange of the 
dilaton, giving rise to a non-linear superposition principle for solutions to Einstein's 
equations. The fact that the chiral primary operators have a ring structure (the 
chiral ring) tells us that is the expected behavior from the J\f = 4 supersymmetric 
field theory. Thus, it should be possible in principle to find fully non-linear solutions 
to the Einstein equations that respect the supersymmetries. 

For large quantum numbers, one enters into the very non-linear regime of the 
classical gravity theory and it is possible to find new topologies and configurations 
that at first sight do not seem to be described by a gravitational background. 

The simplest such configurations arise for BPS states with J ~ A^, such that the 
angular momentum is only happening in an 50(2) subgroup of the 5*0(6) R-charge. 
These new semiclassical configurations correspond to giant gravitons that respect half 
of the supersymmetries. Giant gravitons are D3-brane configurations in AdS^ x 
that respect some of the supersymmetries, in this case the same supersymmetries 
that ordinary gravitons respect. Giant gravitons were introduced in as the 
gravitational dual explanation for the fact that traces oi N x N matrices Z, are not 
algebraically independent. Namely, that 

tr(Z"+i)~5^n^^(^') (3.1) 

i 

with ki < N. We are also required to have Yl^i ~ This means that the 
description of states as a Fock space in terms of traces is truncated. A more complete 
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analysis shows that there are two types of giant gravitons. Some grow into the 5*^ 
direction, and some grow into the AdS direction pO| . 

Their dual description in terms of field theory operators was conjectured in the 
works [H, 21]. The giant gravitons growing in the 5*^ are related to sub determinant 
operators. These sub determinants are the coefficients of the expansion of det{Z + x) . 
These are special cases of Schur functions related to completely antisymmetric repre- 
sentations of U (N). The giant gravitons growing into the AdS were then conjectured 
to be related to Schur functions of completely symmetric representations of U{N). 

As we saw in the previous section, it was convenient to introduce determinants to 
build coherent states. Here we see that determinants are also interesting because of 
their relation to giant gravitons. We will comment more on this relationship further. 

The description of was shown to be equivalent to an integer quantum hall 
droplet picture of two dimensional free fermions in the lowest Landau level with 
a confining quadratic potential [Q. In the fermion description, the two types of 
giants correspond to the holes and particles of the Fermi liquid. The giant gravitons 
corresponding to particles are naturally related to eigenvalues (this was also suggested 
in ||20|).The edge fluctuations are the ordinary gravitons. The particles and holes 
dissolve into edge fluctuations if their energy is of order ^/N. In this regime, the 
corresponding giant gravitons are very small, and all the relevant physics is captured 
by focusing on a plane wave limit geometry. In this case this would be the maximally 
super symmetric plane wave (See also [p^). 

It is natural to expect that when we place many mutually supersymmetric D3- 
branes moving on AdS^ x on top of each other, that we will be able to replace the 
brane stack by a near horizon geometry free of singularities (these are non-dilatonic 
branes [^, and the AdS^ x global geometry of is a typical such example). 
In the half-BPS case, in the fermion language, placing the fermions as near to each 
other as possible, but away from the edge, produces new topologies of the droplet. 
The same is true if we work with many holes on top of each other. It is also natural 
to expect that we can turn on coherent states of gravitational perturbations of a 
system. These coherent states are edge waves of flnite amplitude. Thus one expects 
that any macroscopic shape of droplets corresponds to a gravitational conflguration. 

This possibility of a simple fermionic description in fleld theory motivated Lin, 
Lunin and Maldacena to try to classify all solutions of supergravity that respect 
half of the supersymmetries of AdS^ x and that have the same isometrics as the 
unbroken symmetries that the dual states preserve [Q. 

The metric needs to have the following symmetries: there is an unbroken 5*0 (4) 
because we chose spherically invariant conflgurations on the 5*^ boundary. One has 
an unbroken SO {4) subgroup of the R-charge (the little group of the highest weight 
state), because we only chose perturbations with J being non-zero in an 5*0(2) 
subgroup of 5*0(6). 

There is an additional unbroken translation symmetry because A = J, so there 
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is one extra generator of the superconformal group that annihilates the state and 
commutes with the two 5*0(4) symmetries. 

The natural ansatz for a space with this symmetry is of the form 

ds^ = A{x)dQl + B{x)dQl - C{x){dT + Vi{x)dx'f + gijdx'dx^ (3.2) 

namely, a fibration with two three spheres ^3 and ^3, and a killing vector dr- This 
leaves a three dimensional space of the x variables where all the non-trivial dynamics 
is happening. After imposing the additional requirement of supersymmetry, the met- 
ric simplifies further, and one finds a preferred coordinate system for the remaining 
three variables, that is described by The full metric is then given by 

ds'^ = -h-'\dt + Vidx'f + h^^dy"^ + dx'dx') + ye^d^l + ye'^d^l (3.3) 

where h~'^ = 2ycosh.G, and V can be solved for if G is known. In the conventions 
of WM, one uses the auxiliary function 



^ = ^tanhG (3.4) 
and it was shown that z satisfies a linear partial differential equation 

dAz + ydy(^^^ (3.5) 

There is a potential conical singularity in the metric at y = 0. This is avoided 
if z\y=o takes the value ibl/2. One then finds that h is finite and that the metric 
is regular. We impose these conditions as boundary conditions of z at y = 0. We 
see that the solution is determined by a two coloring of the x^, plane (the regions 
where z = ±1/2, let us call them black and white). This is very similar to the two 
dimensional picture one obtains from the quantum hall droplet. 

The interesting thing to notice here is that y = corresponds to a degeneration 
of one of the two three spheres to zero size. If one matches the supergravity solution 
to the quantum hall picture, the regions that have particles in them correspond to 
the locus where ^^3 (the sphere of AdS^) has vanishing size (the black regions). 

For the vacuum state, AdS^, x S"^, in global coordinates, we have that 

ds^ = - cosh'^ {p)dt^ + dp^ + sinh^ pd^l + d^l (3.6) 

and the locus where the radius of the ^3 shrinks to zero size is exactly the at the 
origin of AdS^, namely p = 0. 

Going back to the quantum hole analogy, if one considers a probe giant graviton 
of the particle type, one adds a small black disc outside the black region. In this 
disc the same that shrinks at the bottom of AdS, shrinks to zero size. Thus we 
can identify the locations of the AdS giant gravitons with the regions where the 5*^ 
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of the conformal boundary shrinks to zero size. For the other giant gravitons (the 
holes in the quantum hall picture), it is the other that degenerates to zero size. 

A more general metric for the 1/8 BPS states would have lower symmetry. Using 
the AdS/CFT dictionary, the states we need in the field theory are constant on the 
three sphere, as we have discussed previously. This symmetry would correspond in 
supergravity to an isometry of the metric. Thus we would still preserve the SO (4) 
isometry, and this requires the metric to have a three sphere fibration. One would 
also have a Killing vector associated to the BPS bound. 

This would give us 

ds"^ = -A^idt + Vidx'f + Bdnl + gijdx'dx^ (3.7) 

This reduces the full gravitational problem to a system of non-linear partial 
differential equations in six dimensions. Such a possibility has been investigated in 



2^, and in a slightly different context in ||26|. One would expect that there are 
regions where this three- sphere fibration degenerates to spheres of zero size. These 
should be the locations where the giant gravitons that grow into AdS are located. 

3.1 Proposal for comparison of matrix models and gravitational physics 

So far, in the gravity side, we have noticed that the regions where the metric is 
interesting correspond to the locus where giant gravitons that grow into AdS are 
located. Fortunately, these giant gravitons can also be identified with the eigenvalues 
of the matrix model itself under the AdS/CFT dictionary |2y, ^, |[. 

Comparing to our previous analytical results of the CFT side in section |], we 
noticed that for semiclassical states the eigenvalues of the matrices form well defined 
submanifolds in the six flat dimensions. It is natural to identify the singular eigen- 
value distributions obtained from the field theory (in the large N limit) with the 
geometric locus where the size of the S"^ fiber of the 1/8 BPS geometry vanishes. 
This S^ survives in the conformal boundary of AdS^, and it represents the S^ on 
which the field theory has been compactified. 

One also expects that there is locally a similar coordinate to y in equation |3.3| , 
such that ?/ = is the degeneration locus. In the case of the ground state, namely 
AdS^ X S"^, the role of the y coordinate is played by the radial direction of AdS, and 
the locus ?/ = is the bottom of the AdS potential well. 

In the eigenvalue picture, the y coordinate can then be thought of as a transverse 
coordinate to the eigenvalue distribution in flat MP . 

For the case of AdS^ x S*^, the six dimensional space corresponds to the radial 
direction of AdS and the 5*^ together. One can map this to the region outside the 
five-sphere distribution of eigenvalues. It is easy to convince oneself that exciting 
one eigenvalue in a BPS manner removes it from the five sphere and places it outside 
the S^ distribution but not inside. 
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Our proposal is that the support of the eigenvalue distribution represents exactly 
the degeneration locus of the three sphere in the full ten dimensional metric. As such, 
one can compare the degeneration locus of exact supergravity solutions with their 
dual description. The tests can be both qualitative, at the level of comparison of 
gross topological features, and they can also be more quantitative. For that, we need 
a proposal on how to extract metric data from an eigenvalue configuration. 

The natural way to do this is to add small (but not massless) strings in a par- 
ticular region. The length of the string serves as a probe of energy. It is also 
natural to add strings in the eigenvalue problem, similar to the work [Q, by turning 
perturbatively the off-diagonal modes of the matrix model. This gives us, to first 
approximation, the induced flat space metric on the submanifold of the eigenvalue 
distribution. 

One can in principle compare the energies of these two different descriptions of 
the same objects, in the degeneration locus. In the gravity setup, one has to take 
into account that time is warped to make the comparison. The comparison of metrics 
seems to give 

Our analysis suggests the following conjecture: for the 1/8 BPS metrics with 
asymptotic AdS^ x 5*^ boundary conditions, the classifications of metrics can be 
reduced to a six dimensional problem that is a base of an S''^ x M fibration. This 
six dimensional (base) space should be diffeomorphic to an open connected subset of 
with a five- dimensional boundary. The region of that is carved out ends on 
the eigenvalue distribution in the dual field theory representation, and contains the 
region around infinity. The fact that the base is diffeomorphic to a subset of is a 
very strong topological constraint. 



4. The Monte Carlo algorithm 

As we described in section H, we are interested in computing the distribution of par- 
ticles for an auxiliary statistical mechanical system of particles in (or in complex 
coordinates, particles in C^) constructed as follows. We start with a wave function 
ijj of the form 

^ = exp(- J2lx^ + Yl - ^' ^) (4-1) 

i i<j 

where P is a holomorphic function, invariant under permutations of the particles 
labeled by i. In particular, the ground state is the case for P = 1. In general, for 
geometric states we expect that logP(X, F, Z) ~ m J2i ^ogg{xi, i/i, Zi), where m is an 
integer and g is a. polynomial in the complex coordinates with complex coefficients. 
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The wave function leads to a probability distribution 



|^P = exp(-5^C 



- m[\og{g{xi, yi, Zi)) - hg{g{xi, yi, Zi))]) + ^ log \xi 



i<j 



(4.3) 



that is equivalent to a Boltzman distribution for a gas at (3 = 1 and with 



H 



f - - m[log(5f(a;i, y^, Zi)) + \og{g{xi, yi, Zi))] - ^ log \xi - Xj 



(4.4) 



The term ^j(x^ — m\[og{g{xi, yi, Zi)) + \og{g{xi,yi,Zi))] can be interpreted as a 
confining external potential, that is asymptotically quadratic. The term that sums 
log \ xi — Xj p can be interpreted as repulsive logarithmic interactions between the 
particles. We are also required to study the system in the thermodynamic limit 
large, where is the number of particles in the system. 

One expects on general grounds that the system will settle to some (thermal) 
equilibrium between the confining potential, the repulsion of the particles and with 
some extra thermal fluctuations. We want to study the (density) distribution of the 
particles in R^, and we also want to measure the fluctuations around the equilibrium 
configurations. 

The equilibrium configuration is a typical configuration of the particles that 
dominates this thermal ensemble. For coarse grained observables (like the typical 
number of particles in some sufficiently large specified region), a typical configura- 
tion is enough to characterize the ensemble accurately. Also, if one has a sufficiently 
symmetric situation, measuring on a single typical configuration can also let us mea- 
sure thermal fluctuations of the ensemble by using symmetry operations to obtain 
more measurements out of a single distribution. 

Our studies for this paper will only involve a special case of the set of distributions 
above, where g{xi, yi, Zi) = Xi, and m is either zero (the ground state) or some other 
positive integer. 

For the case m = 0, we expect that the distribution will have all particles forming 
a round five sphere shell, at a distance of 



with some thermal fluctuations, and in principle 1/A^ corrections to the radius [§]. 
The thermal fluctuations are also a 1/A^ effect. This is because in the large A^ 
counting of t'Hooft, 1/A^ serves as a measure of h, and our probability distribution 
arises from a quantum wave function. The quantum fluctuations of a wave function 
are typically an effect of order h, and they translate to thermal fluctuations of our 
ensemble. In the gravitational theory, 1/A^ effects are quantum gravity corrections 
to classical results. 



(4.5) 



- 19 - 



Thus, in principle, we can measure some quantum gravity effects/corrections by 
looking at these fluctuations in the dual CFT description, that has been simplified 
to this thermal ensemble. 

We use a standard Metrtopolis algorithm for a Markov process that navigates 
between the different configurations of particles in MP. Between successive configu- 
rations A,B,we always accept the B configuration if Hb < Ha, and if Hb > Ha, we 
generate a random number s and compare s to exp(— — Ha))- If s is smaller 
than this quantity, then we accept the B configuration. Otherwise we reject it. 

Our algorithm to compute the thermal ensemble begins by setting up a random 
distribution of particles in six dimensions. Wc use the number of particles to set 
the size of this distribution. Wc disperse the particles randomly on a box of size 2r, 
where r is the expected radius of the distribution, with each coordinate generated 
randomly in this range. 

We generate new trial configurations by moving one particle at a time, — > 
Sxi. The movement 5xi is generated randomly by varying the coordinates in a range 
±5 for each coordinate. After moving one particle, we apply the metropolis criterion, 
and we cycle between the particles in order. After some number of iterations, we make 
S smaller. This is designed to converge faster to a near equilibrium configuration, 
and then we use the smaller value to ensure that reasonable thermal fluctuations are 
generated somewhat accurately. After some prescribed number of iterations / we 
stop the calculations and look at the final configuration and consider it to be typical. 
We usually have half the iterations at one value of 5, and then we make 5 smaller by 
a factor of 5. 

The typical values of N that we use range between 100 and 20000. We found 
that moves with S of order 2 — 5 give good results. We used typically 300 to 2000 
iterations per particle. Experimentally, we were finding convergence of results after 
about 100 — 150 iterations. Our computer code is a C program. We compiled our 
codes with gcc in various Apple computers. We performed all our calculations with 
double precision, and we used the GNU scientific library random number generators 
to insure good statistical quality of the pseudo-random numbers. 

The typical running time oscillated between 5 minutes and 48 hours depending 
on the number of particles and iterations, as well as on the speed of the computer 
processor that was used. The number of floating point operations that we perform 
scale roughly as IN'^, so making N large is computationally expensive. 

4.1 Taking the large N Limit 

We have found that taking N large in the simulation directly is difficult (computa- 
tionally expensive). If the conflgurations are going to be thermodynamic, it makes 
sense to coarse grain the dynamics in order to increase the effective value of N with- 
out increasing the number of computations significantly. The simplest way to do 
that is to assume that in the simulation each particle counts for M particles near 
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their vicinity, so that the effective is NpartM, where M is a multiphcative factor 
associated to this coarse graining. This is taken care of by increasing the strength of 
the repulsive interactions by a factor of M. This is, we are simulating instead 

exp{-(3H) = exp (- J^i^l - H^og{g) + \og{g)]) + ^ Mlog \xi - x,A (4.6) 

\ i=l i<j J 

Roughly, particle i feels the effect of M particles at the location of particle j, but we 
are not averaging the position of particle i over M particles. We need to verify that 
this is a good approximation for the questions we are asking. 

We do this by comparing distributions generated with the same Neff, but with 
different values of M (the coarse graining factor). If we can show that numerically 
the answers are similar, then we can take N^ff to be large by keeping the number 
of particles fixed and varying M. We expect that this would induce errors of order 
where N is the number of particles in the simulation, not N^ff = NM, because 
we are ignoring the M — 1 particles in the vicinity of the i — th particle that i is not 
being affected by. So long as we can show that our errors are of this order, taking 
Neff to be very large is possible and we can compare to the thermodynamic limit. 



5. Results for sphere and Comparison to theory 

In the following tables, we show the numerical results for simulations. In particular 
we are able to show that in the large N limit we approach the expected theoretical 
result from equation |2.19|. It is also evident that the effects of coarse-graining do not 



affect the calculations of the radius of the sphere substantially, and this allows us to 
extrapolate to relatively large values of N^ff = NM by taking M large. 

Here we also present some results of varying the size of the motions of the 
particles in the Monte-Carlo algorithm. There is no substantial difference between 
runs where we take the particles and move them by bigger or smaller steps. 

We can also plot a typical configuration of particles projects into the 1,2 plane, 
for a simulation with = 1000. The projection looks very round, and it is noticeably 
more dense in the center than in the outer edges, exactly as one would expect. 

We can also show a histogram of the radius of each particle. It is evident by 
sight that the approximate radius of each particle is close to the expected value of 
r ~ 22.36. 

Our results comparing the simulation with coarse-graining to the full simulation 
for large numbers of particles are in the table |I|. 

In the table we calculate the radius of the distribution by using 

N 

i - 
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N 


M 


Neff 


'^dist 


rth 


Relative difference 


100 


1 


100 


7.204 ±0.039 


7.07 


±1.9% 


200 


1 


200 


10.094 ±0.034 


10 


±.9% 


100 


2 




10.068 ±0.052 




±0.7% 


400 


1 


400 


14.212 ±0.031 


14.142 


±0.5% 


200 


2 




14.162 ±0.056 




±0.1% 


100 


4 




14.187 ±0.0465 




±0.3% 


800 


1 


800 


20.038 ±0.014 


20 


±0.2% 


400 


2 




20.025 ±0.021 




±0.1% 


200 


4 




20.0165 ±0.023 




±0.1% 


100 


8 




19.928 ±0.041 




-0.4% 


1600 


1 


1600 


28.315 ±0.013 


28.284 


±0.1% 


800 


2 




28.305 ±0.012 




±0.1% 


400 


4 




28.301 ±0.028 




±0.1% 


200 


8 




28.252 ±0.032 




-0.1% 



Table 1: Comparison of coarse graining to fulf simulation 



and we are showing the theoretical value on the right, a/ Neff/2. We also show the 
relative percent difference (r^j^ — rth)/'rth to the theoretical value. The error bars are 
based on samples of 10 configurations per distribution, where we quote the standard 
deviation of a sample around the mean value. Here, it is useful to notice that just like 
in the case of spherical membranes in Matrix theory [^, there are various possible 
definitions of the radius once the statistical fuzzyness of the sphere is taken into 
account. The one we chose seems hke a sensible definition, as it is given by 

r' - ^tr(X2) (5.2) 

the simplest matrix model correlator that one could use. 

Our error bars are the statistical errors of our sample, but they do not include 
an estimate of systematic errors of the code, or show that correlations between dif- 
ferent samples associated to the same distribution are absent. Our purpose in this 
paper has been to establish that this type of simulation is a viable avenue to un- 
derstand geometrical information using a very simplified matrix model analysis, and 
this requires us to take Ngff as large as possible. 

As can be clearly seen, for the most part the effect of coarse graining is consistent 
with a full simulation, and it tends to lower the value of r^ist with respect to the full 
simulation. This is expected, as in our coarse graining procedure we are ignoring 
the self repulsion of the eigenvalues that have been coarse grained into a single one. 
This is a systematic effect that is of order M/N^ff = 1/N and tends to lower the 
size of the distribution (there is a slightly smaller net effective repulsion). We should 
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M 


'^J I 


^ dist 




Relative difference 


1 


400 


14.213 ± 0.031 


14.142 


±0.5% 


2 


800 


20.025 ± 0.021 


20 


±0.1% 


4 


1600 


28.301 ± 0.028 


28.284 


±0.1% 


10 


4000 


44.689 ± 0.255 


44.72 


-0.1% 


20 


8000 


63.178 ±0.027 


89.44 


-0.1% 


50 


20000 


99.890 ±0.012 


100 


-0.1% 


70 


28000 


118.176 ±0.025 


118.32 


-0.1% 


100 


40000 


141.237 ±0.018 


141.42 


-0.1% 


1000 


400000 


446.663 ± 0.027 


447.21 


-0.1% 



Table 2: Coarse grained distribution for large limit. 



also notice that the size of the error bars in the radius distribution are very similar 
between the coarse-grained distribution and the full simulation. Also, the results are 
close to the expected theoretical result. Looking at the simulations with M = 1, we 
notice that as gets larger, the statistical error bars in the radius tend to go down 
in size. For = 1600, the relative error is 0.1%. 

Table ^ shows a simulation with = 400 where we just change the coarse- 
graining multiplicity M through various values. The relative difference is of order 
0.1% ~ 1/1000 for large values of N^ff. This is smaller than 1/400 by a factor of 
2, a systematic M/Ngfj effect. We can estimate this effect by comparing 
to {N — l)M/2, the first is the effective radius one would get for Neff = NM, 
while the second radius is the mean field theory (the particle i feels the repulsion of 
(A^ — 1)M particles and the scaling of the problem. This shows a relative difference 
of order —1/2N ~ —1/800, exactly like we are seeing in the data. Our results are 
therefore consistent with the theoretical large A^ limit when we take this systematic 
effect into account. 

5.1 Density fluctuations and 1/A^ counting 

As we have shown above, our numerical results seem to be converging to the radius 
of the particle distribution that is obtained from the saddle point approximation. It 
is convenient also to test the density fluctuations along the flve-sphere, to test how 
spherically homogeneous it is. One would also like to understand how thin is the 
sphere, as in the saddle point limit we are suppose to be comparing it to a delta 
function distribution at flxed radius. 

To see whether the system has rotational symmetry, one can for example take 
a collection of points generated this way and project them onto the 12 plane. This 
can give us an idea of how spherically symmetric the distributions are. We present 
here an example with A^ = 2000 in flgure |^. 
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Figure 1: Projection onto the 12 plane of a typical configuration of 2000 particles. 

By inspection, the distribution looks fairly circular. We should also remember 
that this is projection of a 5-sphcrc embedded in flat six dimensions. This gives us 
an image of the 5-sphere as a 3-sphere flbration over a disc. The radius of the 
fiber is r| = (r^ — r^2) terms of the (ru) radial coordinate, and r, the radius 
of the distribution in 6 dimensions. The number of points for each sphere element 
scales with the radius of the 3-sphere, so the density of points at the center has to 
be much higher than at the edges. This is observed in the figure. This makes it hard 
to estimate how spherical the object is, because our sight perceives mostly the edges 
of the disc, where the density of coverage is small. 

To do better work, one would want to count the number of particles in cones 
of given angles. However, this can depend very much on how one chooses to slice 
these angular regions and how one triangulates the six sphere into smaller pieces. 
Instead, as is usual in other situations, one would much rather consider the multipole 
expansion of the distribution. Since the particles are all expected to be essentially at 
the same radius, the multipole expansion will give us directly the angular fiuctuations 
of the density into different multipole moments. 

Here, we will just calculate the quadrupole moment squared of the distribu- 
tion. To ensure that we are doing things carefully, we can expand the square of the 
quadrupole moments on all 2-planes. 

The typical quadrupole moment (along the 12 plane lets say) can be written 
as Q ~ 3?e(tr(Z^)) where we choose Z = Xi + iX"^. We want to normalize the 
quadrupole to refiect just the angles and not the radius, so we should use 

Q = 3fJe(tr(Z7r2)) (5.3) 

Because of spherical symmetry, the average quadrupole will be exactly zero. But 
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N 




100 


2.74 


400 


2.56 


800 


3.28 


1200 


2.82 


1600 


3.21 



Table 3: Quadrupole moment table 



to understand the typical fluctuation, we would want to calculate 

Q2 ~ (tr(ZVr2)tr(ZVr2)) (5.4) 

These are typical correlators in matrix models. In the large limit, formal 
developments around perturbative constructions lead to a 1/iV^ expansion into non- 
planar diagrams, so long as one can think of the solution of the matrix model as 
a resummation of a perturbation series. We should notice that the Monte-Carlo 
simulation is non-perturbative in nature, so it is possible to observe different behavior 
than the expected perturbative result in the simulation. 

In our case, we have that the radius r is of order \fN , so that a correlator ii{ZZ) 
scales like iV^, and tr(ZZ)^ would scale like N'^. If our results can be understood in 
this type of framework, we would expect that the numerical values of various objects 
will have a 1/A^ expansion. The tree level value of tr(Z^) is zero, so for |tr(Z^)p 
we expect a result which is suppressed with respect to A^*^ by a factor of A^^. If 
we correct for the normalization of the radius, as in equation |5.4| , so that we can 
compare different droplets of varying A^, we expect that these correlators are all of 
order one. This can be tested in our numerically generated distributions for various 
values of A^. This is shown in table |^. 

The Q is generated by averaging over quadrupoles in various planes. For each 
configuration of eigenvalues, there are 20 linearly independent quadrupole moments. 

In practice we calculate the six numbers 

S,, = tr(X2 - X|)/r2 (5.5) 

and the fifteen numbers 

Q., = tr(2X,X,) = tr[i ((X, + X,f - (X, - X,f)] (5.6) 

and then we average their squares, to obtain Q^. These 21 numbers contain all 
20 independent quadrupole moments of a single distribution of particles. We then 
average these over ten different distributions with the same number of particles. 

If instead we had a gas of uncorrelated particles at constant radius, the total 
quadrupole moment would scale with N"^/"^ where N is the number of particles. The 
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matrix model behavior signals a highly correlated system compared to a free gas. 
In this sense, it can be considered as a type of liquid (the same gas model in two 
dimensions is the particle distribution associated to a quantum hall fluid). 

Indeed, our table |] shows that the total quadrupole moments do not vary sub- 
stantially between various values of the number of particles and that the numerical 
values seem to agree with the expectations of large matrix model statistics. Also 
notice that any growth in the quadrupole moment with increasing is at best 
marginal. We also don't have error bars associated to these numbers Q^: we don't 
have enough information to compute them. The table should be interpreted as pro- 
viding numerical evidence for the expected large A^ scaling of a matrix model. 

It should be noticed that the multipole expansion in angular coordinates would 
give us correlators of the form 

~ (tr(Z"/r")tr(Z"/r")) (5.7) 
while some multipole correlation functions would be of he form 

~ (tr(Z'™/r'")tr(Z"/r")tr(Z7rP)) (5.8) 

and in principle can be compared to three point functions in the conformal field the- 
ory, and via the AdS/CFT, to supergravity. These three point functions for half-BPS 
states were first studied in detail in [Q, and suggested various non-renormalization 



theorems. These have been analyzed in more detail in pl[] . 

To do a useful numerical comparison we expect to need a large statistical sample, 
as these normalized correlators will be expected to be of order 1/A^. Such a statistical 
sample is beyond the scope of the present paper. This is currently under investigation 



Similar to the angular fluctuations of the density described above, one can also 
slice the points radially and try to understand the radial distribution of particles 
in more detail. In particular, one would like to determine if the radial distribution 
of a single particle is approximately Gaussian or not. Naively, non gaussian be- 
havior should be suppressed by some (possibly fractional) powers of 1/A^ and could 
correspond to some non-trivial quantum gravity effect. 

We find that the typical width of the radial distribution of particles is always 
of order 1 (0.86 in the example), while the radius is of order \fN . We also find no 
noticeable deviation from a gaussian distribution for the distribution of the radii for 
a simulation with A^ = 2000. this is encoded in the figures ^ and ^. 

Finally, we can also look at the radii of the individual particles in the order they 
were generated. This is shown in figure ^ 

In this sense, in the large A^ limit the numerical simulation of the density ap- 
proaches a delta function, as the width of the distribution of particles in the radial 
direction goes down relative to the radius of the five-sphere. 
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Figure 2: Binned radial data for one sample of = 2000 on 50 bins between 27 and 37, 
with an expected radius r 32 




Figure 3: Comparison between simulation and a Gaussian fit by bin. 




500 1000 1500 2000 



Figure 4: Typical radii for a configuration of individual particles (the data is by particle 
number) 

6. One hole states 

So far we have done a systematic analysis of the basic properties of the ground state 
wave functions for various A^. Our results so far suggest that it is possible to work 
a detailed comparison between the numerical results of a Monte-Carlo analysis and 
the expectations that one has from a theoretical analysis of a formal expansion. 

We would now like to analyze other configurations that in gravity could corre- 
spond to a non-trivial topology change. In particular, it is interesting to analyze the 
typical distribution of eigenvalues for the simplest LLM geometry whose topology is 
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different than the ground state. 

In particular, one can consider geometries that correspond to an annulus in the 
LLM plane. These are conjectured to be the result of condensing H N maximally 
giant gravitons associated to some particular half BPS orientation. A single such 
giant graviton is believed to be given by the operator det Z, and H such gravitons 
would correspond to the operator det^ Z. In the matrix model of eigenvalues, these 
would correspond to a wave function 

V^H~4det(Z)^ (6.1) 

and to an associated thermal ensemble 

I^P = exp(^(-f2 + H\og{{x'f + (x'f)) + 5^1og(|f. - x,\) (6.2) 

Thus the particles are repelled from the locus = (0,0). In two dimensions, 

such wave functions correspond to an inner circle of the annulus with a radius that 
scales with \/H, independent of A^. 

Also, from the LLM intuition, one expects that in the LLM plane the region 
with eigenvalues (particles) is the degeneration locus of the 5*^ of the boundary. At 
the edges of the droplets, both the of the boundary and the associated to 
being written as a fibration over a disc vanish. Thus we expect that if our particle 
distributions are to match the LLM topologies, then the projection of the particles 
to the (1, 2) plane will form an annulus, and at the edges, the size of the sphere in 
the orthogonal directions should vanish. As described previously, matching to an 
exact LLM metric is difficult, as it is not clear what is the coordinate change that 
relates the eigenvalue distributions to the LLM coordinates, and this map could be 
rather complicated. However, matching the topologies of the degeneration locus of 
the boundary should be straightforward. 

As can be seen from figure ||, the distribution of particles projected on the plane 
forms an annulus shaped distribution of particles, exactly as expected from a naive 
analysis of the topology of such an LLM droplet. In this sense, we are observing 
directly a topology change. However, the simulations do not show the expected 
closing of the distribution of particles into a donut, but rather there appears to be 
a boundary to the surface geometry that the particles describe. We do not know 
at this moment if this is a feature of the physical system we are studying, or if we 
need to improve the code to focus on the region where the deviation from our naive 
intuition is taking place. This is currently being investigated. For an alternative 
picture, we can plot the radius on the 12 plane with respect to the radius in the 3456 
plane, as shown in figure ^ 

If it is a new feature of the system, we need to explain the physics associated to 
it, and we may need to rethink some issues regarding what we mean by geometry 
associated to these droplets. We believe this is a very interesting result of the simu- 
lations we have done. Below, we will present evidence that this type of simulation is 
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Figure 5: The plot shows the average radius of the particles inthe 3456 directions for 
binned data in the 12 plane. The data shows a simulation with 2000 particles and H = 300, 
in a 40 X 40 grid. 




Figure 6: ri2 vs ri23456 for an = 2000 particle distribution, with H = 30. Notice the 
low number of points for ri2 large, relative to ri2 small: this is an effect of phase space for 
each type of configuration. 

behaving like expected in terms of the AdS intuition for the size of the hole in the 
center, giving us confidence that the distribution ending "in the wrong place" may 
be a physical feature that has to be considered seriously. 

It is customary to compare LLM droplets directly with distributions of fermions 
associated to the lowest Landau level in two dimensions for a given wave function. 
That is the simplest description of all half-BPS states |l4l. Similarly, various 



studies of superstar geometries exploit this correspondence ||3^ . 

For such wave functions, the size of the center hole in the Fermi-liquid picture is 
determined by the number of hole-particles that one places at the center of the ring. 
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N 


R 


400 


8.45 ± 1 


600 


9.26 ± 1 


800 


10.0 ± 1 


1200 


11.1 ± 1 


1800 


12.57 ± 1 


2400 


13.34 ±1 



Table 4: The above shows simulations with varying numbers of particles and fixed number 
of maximal giant gravitons {H = 30). The right hand side is computed by averaging the 
five lowest values of i2i2, the radius in the 12 plane. 



and this size is independent of the droplet. 

In our numerical simulations, we see a dependence of the size of the hole in terms 
of the total number of particles. This is shown in the table |. The way we determine 
the inner radius of the distribution is by averaging the five smallest values of R12. 
The choice of 5 is arbitrary and it is done to reduce large fluctuations of a single 
particle. The error in the simulations seems to be dominated by systematic errors, 
and we don't know how to estimate them, but we believe they can be bound by 
±1 on all entries. Changing the coarse graining factor has an effect on the values 
obtained, and it is not clear how to compare different distributions associated to the 
same number of particles A^^e//, unless we put large systematic error bars of that 
size. Also, we have no good theoretical understanding of the fluctuations of the 
distributions. We can not compare results to a large N limit saddle point either, as 
we have no theoretical prediction of what the saddle point will look like. 

We can analyze the data In a log- log plot, as shown in figure]^ to test scaling, 
with N at H fixed. The data can be fit linearly to very good accuracy to a slope of 
0.26, which seems very close to 1/4. 



★ 



G;25 6.5 6.75 7 7.25 7.5 7.75 



Figure 7: Log-Log data of size of central hole, as read from table 4. 

Scaling of the saddle point suggests that the hole radius should be fixed relative 
to the size of the distribution if we keep H/N fixed. In this fashion, we expect that 
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'"in ~ g{H/N) is a function of H/N. For iJ/A^ small, we seem to obtain simple 
scaling behavior. 

Thus, we expect that the radius of the hole is roughly given by r ~ {HNY^'^, by 
guessing a simple rational exponent near the slope we found. 

This result does not mean that the comparison with free fermion droplets are 
wrong. Rather, it shows that the LLM plane can be embedded very non-trivially 
into the M^, with an as- yet unknown change of variables determining this embedding. 
It is also worth pointing out that the LLM metric is not associated to a flat metric 
on the LLM plane, but that there is also a warping by a function h"^, as written 
in equation |3.3| . As discussed in section ^, matching to an LLM geometry needs to 
resolve these changes of coordinates. 

This factor of N'^^^ is reminiscent of the scaling of the string length. Indeed, 
we can motivate this comparison in a straightforward fashion. If we assume that 
we have H D3-branes, where H << N, we can think of the H branes locally as 
flat D3-branes embedded in fiat space, with small back reaction far away from the 
branes. We would find then that the near geometry will start differing from a fiat 
metric in a region of size H^^^. 

This is because for H D3-branes, the metric of the geometry of the D3-branes 
in fiat space is given by p4 



ds' ~ r'^'dxl + f'\dr^ + r^dnl) (6.3) 

with 

+ (6.4) 

where r is a radial direction transverse to the D-brane. The region where / differs 
substantially from one is of order r ~ H^/^^ in string units. To be more precise, 
the region where we expect a discrepancy from fiat space is of order H^/^ in Planck 
units. Indeed, when we express a' in terms of Ipi, the factors of g'^ will cancel. This 
can be traced to the fact that the tension of the D3 brane does not depend on gg in 
units of Ipi. 

Now, since in the dual AdS geometry to the A/" = 4 SYM with gauge group 
U{N) we find that Ipi scales like N^^'^, we find that having an inner radius that scales 
with (HNy^^ is exactly right to match the intuition from AdS. We find this way 
that the matrix model we have written knows about the Planck scale in a way that 
can be measured geometrically. 

This gives us evidence that the numerical results are capturing geometrical data 
associated to such a collection of H D3-branes rather precisely, even if the comparison 
to a particular LLM metric is hard. 

7. Conclusion and prospects 

In this paper we have explored how certain aspects of the AdS /C FT correspondence 
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can be studied numerically using Monte Carlo methods. We believe that we have 
shown that this is feasible: we can explore the large limit effectively with modest 
computational resources. The main simplification that we have used is that at strong 
coupling and for certain states only constant matrix fields contribute to the dynamics, 
and moreover the dynamics reduces to configurations of commuting matrices. Thus, 
for large matrices of rank N, the number of degrees of freedom that are relevant are of 
order A^, as opposed to A^^. We also have a list of exact wave functions to simulate. 
Thus, we do not need to study the time evolution of degrees of freedom, but 
instead we can concentrate on thermodynamic properties of the degrees of freedom. 

In particular, the most important description of the system under study is in 
terms of a distribution of particles in six dimensions and the geometric properties of 
these distributions can be correlated with certain geometrical degenerations of the 
higher dimensional AdS dual geometry. 

To the extent that there are exact analytical results for these distribution, the 
numerical data seems to match the large A^ limit very precisely. This is possible only 
for the ground state, where we have maximal rotational symmetry and is the only 
known distribution of particles in the large A^ limit. 

The distribution in the saddle point limit is a delta function. Since the radial 
fluctuations of the distribution are of order one, while the radius of the distribution is 
of order A^^/^, the distribution approaches a delta function in the limit. It is expected 
that there might be 1/A^ corrections on the radius depending on what observables 
are chosen to deflne the radius itself. Many of these agree at large A^, but differ by 
various 1/A^ corrections. 

It is natural to imagine that these differences might account for various quantum 
gravity effects, as different probes will be sensitive to different deflnitions of the 
radius. 

In some sense, the simulations we are doing encode all the quantum gravity 
corrections to the geometry. However, the dictionary between our simulations and 
gravity needs to be improved. 

We flnd that it is possible to consider other spacetime topologies. Guided by 
the intuition of the LLM metrics [1^ and the construction of giant graviton state 
operators in the CFT [|r^, we have been able to consider wave functions that give 
rise to measurable topology changes in the distributions of particles. In particular, 
we were able to show with these conflgurations that the matrix model of six matrices 
is aware of the Planck scale size of the flve sphere of the dual AdS^ x geometry, 
and therefore we have a measure of a relevant non-trivial scale in the problem. 

It is also worthwhile to point out that the size of the sphere in fleld theory units 
is of order V^, and that each particle occupies a typical volume of size A^^/^ ~ 
(A^-'^/4)^A^^/^ ~ N^^^lpi, much larger than a Planck scale volume. This shows that 
the AdS/CFT correspondence sees locality from a very different perspective than 
other discretizations of gravity, where the Planck scale is the scale of granularity of 
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the discretization [|T|. Instead, we are able to see below the scale of discretization 
by studying collective effects and the Planck scale is related to studying the back 
reaction of the geometry to a heavy membrane. As a curiosity related to this problem, 
in figure |^ the particles accumulate densely in the region where geometry is behaving 
differently and it suggests that the density of resolution compensates for the scale of 
interesting features below the typical scale of granularity. 

We think it is important to emphasize that we are also very puzzled by the 
fact that the particle distributions are not behaving like we expected from a simple 
comparison to an LLM metric intuition. This discrepancy might be due to subtle 
issues regarding the simulation itself. On the other hand it is possible that we have 
not identified the precise correspondence to geometry in detail. Deciding between 
these options deserves further study. 

We believe it is possible to make various measurements of correlators numerically 
and to compare them with the non-renormalization theorems for three point functions 
of BPS states. This will require studying large samples of statistically independent 
configurations, necessitating an improved study of error analysis to complement the 
relatively crude estimates presented here. 

It is straightforward to generalize our calculations to orbifold setups, following 



1^ and to more complicated LLM geometries. In particular, one can imagine that 



once it is better understood how to correlate the properties of the distributions to the 
metric aspects of the LLM geometries, it will be possible to study directly topology 
changing transitions, like the ones considered in [p^ . 

Another important aspect that might be studied with these methods is to con- 
sider a probe giant graviton that grows into the AdS region |T^, Studying the 
corresponding wave functions in detail could give us a precise coordinate map be- 
tween the CFT field variables and the radial direction of the AdS geometry. 

One can also imagine using these numerical techniques to explore field theories 
that are not well understood, like the CFT duals associated to AdS^xS'^ or AcIStxS^. 
A recent proposal states that BPS states for the field theory duals are also described 
by a system of commuting matrices |^ . Unlike the case of A/" = 4 SYM, there is no 
measure that one can derive from first principles to study these objects as we have 
described in this paper. Numerical techniques might be very useful in making checks 
of a guess for that measure. 

Finally, it would extremely interesting if one could extend these techniques fur- 
ther to study string states directly. This involves introducing the off-diagonal modes 
systematically into this setup, as well as the other spherical harmonics of the quan- 
tum fields on the sphere 0. Some partial success has been achieved for adding strings 
to the geometric ground state in [Q. However, one can study the same problem for 
other geometries where analytic string solutions are not available. 

To conclude, we believe that the prospects for these types of simulations we have 
presented in this paper is full of open questions that can be addressed numerically 
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with modest computer capabilities. It is our hope that numerical simulations might 
serve us a source of intuition for effects that we can not calculate otherwise. 
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A. Holomorphic wave functions 

We want to consider the following Hamiltonian for particles in 2d dimensions 



where 



fi" = \[\x^-x,\^ . (A.2) 

and wave functions of the form 

^ = MX)P{z) (A.3) 

where we follow the notation of section In this notation ipQ is a gaussian eigen- 
function of the Hamiltonian. 

7/>o~exp(- 5^x^/2) (A.4) 

i 

P{z) is a homogeneous holomorphic polynomial on 3A^ complex variables 

~ xf-^ + ^xf. (A.5) 

The polynomial is invariant under permutation of the Xj. The z are holomorphic 
coordinates of the particles, after we have made a choice of complex structure on 
M^"^. The Hamiltonian is the leading semiclassical approximation of the effective 
Hamiltonian that encodes all BPS states associated to the chiral ring of A/" = 4 SYM 
theory 0. We are interested in studying how good a description this Hamiltonian 
provides for the conjectured states that we have written in section ^. As we have 
seen, we have a ground state wave function ipQ that is gaussian. We want to show 



that the possible corrections to |A.3| , or to |A.1| are very small, and can be ignored. 
When we calculate ^ 

—V. ■ ^'Vi^oP (A.6) 
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there are terms that contain no derivatives of P, terms with only one derivative 
acting on P, and terms that contain two derivatives of P. Let us separate these 
terms. 

For the terms that do not depend on derivatives of P, P is acting just hke a 
constant, and we can use the fact that ipo is an eigenfunction of H to show that these 
terms just give Eqi^q, the energy of the ground state. Now, the terms that contain 
two derivatives of P are necessarily proportional to 

(Vi)2p = (A.7) 

and these vanish because P is holomorphic. Indeed, is is possible to write 

s 

So we are only left to understand the terms that contain one derivative of P. 
These are given by 

^ ( ViP) • (V^o Vi/x^ + /.22 ViV^o) (A.9) 

The second term gives us ipoXi • ViP, and when we sum over i, we get the Euler vector 
for the scaling xj Xxj. Since P is a homogeneous polynomial of degree m, we find 
that P is an eigenfunction of the corresponding scaling operator, with eigenvalue m. 

The only obstruction for this wave-function to be an eigenfunction of the Hamil- 
tonian is to check if 

EViP-Vi/.2 (A. 10) 



vanishes or not. Now, we can easily show that 



so that we can rewrite the expression in such a way that 



I ilj j JU J 



and there is no obvious mechanism for cancellations (except when d — 1). 

What we want to do now is to think of this wave function as a variational 
approximation to an eigenfunction of the Hamiltonan. Seeing as different degree 
polynomials give orthogonal states (they have different 5*0(6) quantum numbers), 
this is a well defined procedure to estimate the minimum energy solution of the 
Hamiltonian with given R-charge quantum numbers. 

What we notice is that these terms that do not cancel can only be important 
when I X ^ X j I IS small. For this situation we can ignore all other variables and 
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concentrate on two of them {a and /? lets say), and expand P in taylor series about 
their middle point z^. 
This is, we define 

= + Sz,z^ = -Sz (A. 13) 

If we define the following function of and Sz with all other particles at fixed 
position, f{z^, Sz) = P(z", zf^), then / ~ P{z^ + 6z, z° - 6z) = P{z^ - 6z, z^ + 6z) 
because of the symmetry properties of P. From here we get f{z^, Sz) = f{z^, —6z). 
Also, from the change of variables we find that Vq — V/3 ~ V^^. We want to 
understand / near the region where 5z — 0. Prom the symmetry properties, the first 
derivative of / vanishes at Sz — 0, and we need to expand / (the function P) to 
second order in Sz. We find then that the possible singular correction to the energy 
is finite (of order f" / f). 

In this sense, the term that we can not control easily is bounded, and can 
be considered as a small perturbation near the region when the particles coincide. 
Elsewhere, the wave function is an approximate solution of the time independent 
Schrodinger equation. For this last argument we have not exploited any additional 
properties of P other than the symmetry under exchange of particles. 

What we want to do now, is to be more careful about the second order expansion 
of what we have labeled /. Indeed, a Taylor expansion of / shows us that 

f{z\ Sz) ~ /(^o, 0) + ^Sz'Sz^f,ij{z^, 0) + . . . (A.14) 

where the partial derivatives in / are taken with respect to the Sz variables. 
We get then that 

^ = /,,(^o,0)5^^- (A.15) 

and that we need to evaluate the average contribution to the energy by integrating 
this expression in the relevant region where this term is important, namely 

/ (A.16) 

We notice that the integral over Sz an. a. small symmetric disc vanishes unless i — j 
in the expression above, and that spherical symmetry guarantees that this is pro- 
portional to fii{zQ, 0) ~ '^Izf- ^^'^ we go back to our expression for P, and we 
remember that P is holomorphic, so that / is also a holomorphic function of Sz once 
we have chosen a complex structure for the This means that for P holomorphic, 
the term that is hard to compute averages to zero in the relevant region where it is 
important for the case of holomorphic functions. 

This means that the energy of the state has to be given by 

E-^EQ + m (A.17) 
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where m is the degree of P and E is extremely close to Eq + m. 

In the fully supersymmetric problem the corresponding BPS states are supposed 
to have energy given exactly by E = Eq + m, and all other states with the same R 
charge have greater energy of order 1. This implies that corrections from ignoring 
supersymmetry are small, and that we can assume that in the worst case scenario, 
that the wave functions are not exact solutions of the supersymmetric problem, then 
the wave functions we have written are very good variational approximations to those 
wave functions. 
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